ELSEVIER 


Journal of Power Sources 164 (2007) 730-741 


Available online at www.sciencedirect.com 


ScienceDirect 


JOURNAL OF 


www.elsevier.com /locate /jpowsour 


Numerical prediction of the effect of catalyst layer Nafion loading 
on the performance of PEM fuel cells 


Chin-Hsiang Cheng**, Hung-Hsiang Lin®, Guang-Jer Lai? 


à Department of Aeronautics and Astronautics, National Cheng Kung University, No. 1, Ta Shieh Road, Tainan 70101, Taiwan, ROC 
> Department of Mechanical Engineering, Tatung University, Taipei 10451, Taiwan, ROC 


Received 20 October 2006; received in revised form 19 November 2006; accepted 20 November 2006 
Available online 3 January 2007 


Abstract 


Numerical simulations are performed to investigate the effect on PEM fuel cell performance of Nafion loading in the catalyst layer. The 
investigation also considers variations of geometric parameters. A model that accounts for the volume fractions of Nafion, the solid catalyst 
particles, and the void space inside the catalyst layers is incorporated into a three-dimensional computational fluid dynamics code, capable 
of resolving three-dimensional mass, momentum, and species transport phenomena as well as the electron- and proton-transfer processes in a 
PEMFC. Numerical results are first compared with experiments, showing close agreement between predictions and measurements. A parametric 
study on the effects of Nafion loading and geometric parameters variation is carried out to evaluate the performance of PEMFC for various parameter 


combinations. 
© 2006 Elsevier B.V. All rights reserved. 
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1. Introduction 


Fuel cells are an important technology for a potentially wide 
variety of applications including auxiliary power, transportation 
power, and stationary power for buildings and other distributed 
generation applications. These applications are encountered in 
a large number of industries worldwide. 

In the electrochemical reactions of the proton exchange mem- 
brane fuel cells (PEMFC), hydrogen at the anode provides 
protons, freeing electrons in the process that must pass through 
an external circuit to reach the cathode. The protons, which 
remain solvated with a certain number of water molecules, 
migrate through the membrane to the cathode to react with oxy- 
gen and the returning electrons. Water is subsequently produced 
at the cathode. In a PEMFC, Nafion is used as the membrane 
placed between the anode and the cathode. Furthermore, it is 
also employed as an important component of the active layer in 
the catalyst layer. The content of Nafion in the catalyst layers 
affects simultaneously the gas permeability, the catalytic activ- 
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ity and the ionic resistance. Hence, an optimal Nafion content 
in the catalyst layer is definitely desired for good performance. 

Introducing Nafion into the catalyst layer makes the layer 
active in three dimensions because the protons can then be trans- 
ported throughout the entire layer. Boyer et al. [1] reported that 
the ionic conductivity of the catalyst layer increases with the 
volume fraction of Nafion in the composite mixture, and the 
authors derived several expressions for determining the opti- 
mum thickness, Pt loading and Nafion content for the electrodes 
by considering these variables as the single variable in these 
expressions. According to the reports given by Passalacqua et 
al. [2,3] and Antolini et al. [4], it is found that since Nafion is 
an electronic insulator, it decreases the electronic conductivity 
of the catalyst layer, and may cause some catalyst particles to be 
electronically isolated. These isolated particles will not be able to 
participate in any electrochemical reactions. In addition, Nafion 
will change the porosity, permeability, and hydrophobicity of 
the catalyst layer. Therefore, the Nafion content in the catalyst 
layer needs to be precisely controlled in order to achieve the best 
balance among all these influencing factors. Qi and Kaufman [5] 
simply mixed the carbon-supported catalysts with Nafion with- 
out adding any extra material. In this report, it was observed that 
the performance increased with Nafion content up to 30%, and it 
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water activity 

mass diffusivity (m? s7!) 

Faraday constant, 96487 C mol`! 

height of gas channel (m) 

local electric current density (A m7?) 
local ionic current density (A m~?) 
electric current density (A m~?) 

transfer current density (A m7?) 
permeability (m7) 

gas channel width (m) 

gas channel length (m) 

width of a module (m) 

molecular weight (kg kmol~!) 

coordinate in the direction normal to the surface 
pressure (Pa) 

power density (W m~?) 

concentration parameter 

universal gas constant, 8.314 J mol`! K7! 
source term in species equation 

thickness of catalyst layer (m) 

thickness of gas diffusion layer (m) 
thickness of membrane (m) 

thickness of carbon plate (m) 
temperature (K) 

gas velocity vector (ms~!) 

volume occupied by Nafion in catalyst layer 
volume occupied by solid catalyst particles in cat- 
alyst layer 

total volume of catalyst layer 

volume of void space in catalyst layer 
inlet gas velocity (m s7!) 

cell voltage (V) 

mole fraction 

Cartesian coordinates 

mass fraction 


Greek letters 


transfer coefficient for the reaction 

Bruggemann coefficient 

porosity of gas diffusion layer 

volume fraction of Nafion in catalyst layer 
volume fraction of solid catalyst particles in cat- 
alyst layer 

volume fraction of void space (porosity) in cata- 
lyst layer 

electric potential (V) 

phase potential (V) 

ionic conductivity of membrane (Q7! m7!) 
ionic conductivity of Nafion (Q7! m7!) 
membrane water content, 

kmol H20 (kmol $O3~)~! 

gas channel width fraction 

density of oxygen (kg m~?) 


OGpL electronic conductivity of gas diffusion layer 
(Q7! m’) 

Op electronic conductivity of carbon plate 
(Q7! m7!) 

Oscat electronic conductivity of solid catalyst particle 
(Q7! m7!) 

T fluid stress (N m7?) 

la source term in electronic conduction equation 

Subscripts 

a anode 

c cathode/gas channel 

Cat catalyst layer 

eff effective 

GDL gas diffusion layer 

i inlet 

m membrane 

N Nafion 

o outlet 

P carbon plate 

ref reference 

S solid catalyst particles 


decreased afterwards. Their experimental data showed that the 
best performance was achieved with Pt loading of 0.2 + 0.05 and 
0.35 + 0.05 mg cm~?, for 20 and 40% Pt/C support, respectively. 
In addition, Lufrano et al. [6] investigated the effects of Nafion 
loading in electrodes on the performance of supercapacitors by 
electrochemical impedance spectroscopy. The capacitors, with 
lower Nafion content but having higher conductivities and better 
capacitance performances, were obtained. 

Numerical methods have been employed to study the physi- 
cal phenomenon of interest in the fuel cells since these methods 
require relatively lower cost for analysis than the experimental 
ones. Furthermore, they in general provide a sufficient amount of 
data not easily obtainable with measurements. Thus, important 
natures of the flow and the thermal fields might then be observed 
and discussed extensively by numerical studies. To name a few, 
Wang and Savinell [7] performed a simulation of the effects of 
porosity and thickness of the catalyst layer, Pt loading, and CO 
poisoning on fuel cells. Yi and Nguyen [8] used the numeri- 
cal methods to solve a two-dimensional single-phase PEMFC 
model with interdigitated flow channels so as to evaluate the 
effects of inlet and exit pressures, gas diffusion layer thickness, 
and carbon plate width on the performance of PEMFC. Gurau 
et al. [9] developed a half-cell model to deal with the transport 
phenomena at the cathode of the PEMFC. The research explored 
the effects of the porosity of the gas diffusion layer on the perfor- 
mance of fuel cells. More recently, Lin et al. [10] investigated the 
optimization of the proton exchange membrane fuel cell. A sim- 
plified conjugate-gradient method (SCGM) is employed to seek 
the optimal combination of the design parameters, including the 
channel width ratio, the porosity of GDL, and the porosity of 
the catalyst layer. 
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As to the three-dimensional models, Hontanon et al. [11] 
adopted a commercial computational fluid dynamics code (FLU- 
ENT) to predict the effects of a group of geometrical parameters 
on the polarization performance of the fuel cells. Berning and 
Djilali [12] used a three-dimensional physical model to evaluate 
the performance of the fuel cells at different operating condi- 
tions. Their research showed that an increase in the fuel cell 
temperature could raise the cell exchange current density, ionic 
conductivity of proton exchange membrane, gas diffusivity, and 
the operating voltage and then lead to an increase in the cell 
performance. 

In spite of the valuable experiments performed by many 
groups [1—6] to investigate optimal Nafion content in the catalyst 
layers, theoretical information on the effects of Nafion loading is 
relatively scarce. Reliable information about the fuel cell behav- 
ior may be collected by experimental or numerical methods. 
However, unfortunately, so far only a limited number of theoret- 
ical reports have been presented. Among the very few theoretical 
works, Song et al. [13] proposed an optimization problem for 
the cathode catalyst layer performance by assuming that the 
Nafion content and Pt loading are uniformly distributed in the 
cathode catalyst layer. From this, the optimal Nafion content, Pt 
loading, and layer thickness were determined. In a recent report 
[14], the same group of authors studied the optimal distribu- 
tions of Nafion content and Pt loading across the entire thickness 
of cathode catalyst layers, which was obtained by maximizing 
the cell current density at a given potential. When optimization 
was performed with either the Nafion content or Pt loading as 
the single variable, authors found that the optimal distribution 
is a linearly increasing function from the gas diffusion layer 
side to the membrane side. Secanell et al. [15] addressed the 
same topic using a two-dimensional model and solved for the 
catalyst layer composition using a multivariable optimization 
technique. 

With the help of theoretical study, the basics of the improve- 
ment of fuel cell performance by Nafion loading might be 
clarified and the optimal Nafion content could be determined. In 
these circumstances, the present study is aimed to investigate the 
influence of Nafion loading in the catalyst layers on the behavior 
of a proton exchange membrane fuel cell by numerical simula- 
tion. In this study, a self-developed subroutine, which is used to 
deal with the volume fractions of Nafion loaded, solid catalyst 
particles, and void space inside the catalyst layers, is incorpo- 
rated with a three-dimensional computational fluid dynamics 
code, which is capable of dealing with the three-dimensional 
mass, momentum, and species transport phenomena as well as 
the electron- and proton-transfer process taking place in the 
PEMFC. 

Numerical predictions of the Nafion loading effects by using 
a full-cell three-dimensional model are presented. In addition, a 
22 cm x 22 cm single-cell fuel cell module with 14cm x 14cm 
active area has been installed and tested. Numerical results are 
compared with the experiments for validation. A comprehensive 
parametric study of the effects of the parameters including vol- 
ume fraction of Nafion in the catalyst layer, volume fraction of 
void space in the catalyst layer, gas channel width fraction, thick- 
ness and porosity of gas diffusion layer, and membrane thickness 


Carbon plate 
RJ Gas diffusion layer 
E=] Catalyst layer 
RS PEM membrane 


(a) 


© arse 


Fig. 1. Schematic ofa single-cell proton exchange membrane fuel cell. (a) Single 
cell with parallel straight channels and (b) flow channel geometry. 


has been further carried out to evaluate the performance of 
PEMFC under various combinations of these parameters. 
Schematic of a single-cell PEMFC with parallel straight 
channels on both the anode and the cathode carbon plates is 
shown in Fig. 1. The single-cell PEMFC consists of a carbon 
plate, a gas diffusion layer (GDL), a catalyst layer, for each of 
the anode and the cathode, as well as a PEM membrane at the 
center. Here the physical and geometrical parameters of the cell 
shown in this figure may be divided into two groups: the fixed 
and the varied parameters. For these fixed parameters, their val- 
ues are kept constant throughout the study. On the other hand, 
as mentioned earlier, the volume fraction of Nafion in the cat- 
alyst layer (€N,cat), the volume fraction of void space in the 
catalyst layer (€v,cat), the gas channel width (lc), the thickness 
and porosity of the gas diffusion layer (tGp_ and eGpL), and 
the membrane thickness (tm) are regarded as major factors and 
are varied in this study in order to evaluate their effects. These 
parameters are referred to as the varied parameters. Values of the 
fixed parameters are listed in Table 1, and the variation ranges 
of these varied parameters are given in Table 2. Also provided in 
Table 2 are the values of the varied parameters for the base case, 
which is defined to be the test problem for the parametric study. 
Note that in this study, the operating temperature and pressure 
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Table 1 
Fixed parameters in numerical simulation 
Parameter Symbol Value 
Operating temperature (K) T 323 
Inlet and outlet pressure (Pa) Pi, Po 101,325 
Gas channel length (m) lL 5.0 x 107? 
Width of a module (m) lp 1.0 x 1073 
Thickness of catalyst layer (m) tcat 1.0 x 1075 
Thickness of carbon plate (m) tp 2.0 x 1073 
Height of gas channel (m) h 1.0x 1073 
Bruggemann coefficient of gas XGDL 1.5 
diffusion layer 
Bruggemann coefficient of catalyst XCat 1.5 
layer 
GDL permeability (m?) kop 1.76 x 107!! 
Catalyst layer permeability (m?) Keat 1.76 x 1071! 
Anode inlet gas velocity (m s7!) Va 0.3 
Cathode inlet gas velocity (m s7!) Ve 0.5 
Anode inlet mass fraction of H2 Yp, i 0.445 
Anode inlet mass fraction of H20 Ym,o,i 0.555 
Cathode inlet mass fraction of O2 Yo,,i 0.212 
Cathode inlet mass fraction of N2 YN, i 0.709 
Cathode inlet mass fraction of H2O Ymo,i 0.079 
Operating voltage (V) v 0.7 
Anodic charge transfer coefficients Qaa 0.5 
for anode reaction 
Cathodic charge transfer coefficients Qac 0.5 
for anode reaction 
Anodic charge transfer coefficients Qca 1.5 
for cathode reaction 
Cathodic charge transfer coefficients Ace 1.5 
for cathode reaction 
Reference transfer current density at Jaref 9.0 x 108 
anode (A m~?) 
Reference transfer current density at Jeet 250 
cathode (A m7?) 
Anode concentration parameter Ta 1.0 
Cathode concentration parameter Te 1.0 
Electronic conductivity of gas OGDL 300 
diffusion layer (Q7! m7!) 
Electronic conductivity of carbon op 4000 
plate (Q27! m7!) 
Electronic conductivity of catalyst OS,Cat 300 
layer (Q7! m7!) 
Ionic conductivity of Nafion TYN,Cat 25.56 
(27! m7) 
Table 2 
Varied parameters in numerical simulation 
Parameter Symbol Base case Variation range 
Volume fraction of EN,Cat 0.3 0.1-0.8 
Nafion in catalyst 
layer 
Volume fraction of void = evcat 0.112 0.05-0.65 
space in catalyst layer 
Volume fraction of solid — €s cat 0.588 0.05-0.788 
catalyst particles in 
catalyst layer 
Gas channel width (m) lc 5.0x 1074 1.0 x 1074 to 9.0 x 1074 
Thickness of gas tGDL 3.0 x 1074 0.5 x 1074 to 5.0 x 1074 
diffusion layer (m) 
Porosity of gas diffusion eGpL 0.5 0.1-0.9 
layer 
Thickness of membrane fm 1.78 x 1074 0.5x 1074 to 5.0 x 1074 


m) 
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are assigned to be 323 K and 101,325 Pa, respectively. The typ- 
ical value of the operating voltage of the fuel cell is 0.7 V in the 
computation. 


2. Theoretical analysis and numerical methods 


The following assumptions are made in prior to the derivation 
of the governing equations: 


the fuel cell operates at steady state; 

. the gas flows in the gas channels and the porous layers are 
laminar and compressible, and effects of buoyancy on the 
flow motion are negligible; 

3. thermodynamic and electrochemical properties of the gases 
and the solid materials of the fuel cell components are 
assumed constant; 

4. the operating temperature is uniform and constant in the 
solution domain; 

5. the gas diffusion layers and the catalyst layers are homoge- 
neous, isotropic porous media; 

6. owing to high ohmic resistance of the membrane, the inter- 
face between the catalyst layer and the membrane is assumed 
to be insulated against electrons; 

7. water exists in the entire fuel cell only in vapor phase, and 

no phase change process is taken into consideration. 


N. 


2.1. Theoretical model 


The present three-dimensional model involves (1) mass, 
momentum, and species transport phenomena in the gas chan- 
nels and the porous layers, including the gas diffusion and the 
catalyst layers, (2) electrochemical reactions in the catalyst lay- 
ers, (3) electronic conduction in the porous layers and the carbon 
plates, and (4) ionic conduction taking place in the catalyst layers 
and the proton exchange membrane. In addition, the theoretical 
model is applied to evaluating the influence of Nafion loading at 
various combinations of the volume fractions of Nafion loaded, 
solid catalyst particles, and void space in the catalyst layers. 
The governing equations of the theoretical model are described 
briefly in the following. 


2.1.1. Mass, momentum, and species transport (gas 
channels, GDLs, and catalyst layers) 

The governing equations for the gas flows in the gas channels 
and the porous layers are the conservation equations of mass, 
momentum, and species. Application of the flow channels and 
the gas diffusion layers in the fuel cells is to allow the reactant 
gases to diffuse uniformly into the catalyst layers but prevent 
excessive increase in electrical resistance against the released 
electrons from the catalyst layers. The porous layers include 
the gas diffusion layers and the catalyst layers, in which the 
mass, momentum, and species equations for the flows are derived 
based on the non-Darcy law. These equations are given below: 


Mass equation: 


V -(epU) =0 (1) 
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Momentum equation: 


T 
POE Ù 

V. (eoUU) = -eV P + V - (et) + oe (2) 

Species equation: 

V - (epUY;) = V - (D$ËVY;) + Si (3) 


where € and k are the porosity and the permeability of the porous 
layers, respectively; t represents the fluid stress tensor. Fuel gas 
at the anode is considered to contain Hz and H20, and oxidant 
gas at the cathode contains O2, N2, and H20. Y; denotes the 
species mass fraction of species i. Note that for the pure fluid 
flows in the gas channels, the values of £ and k are given with 
£=] and k=. 

Calculation for the effective mass diffusivity of species i 
(D) is based on Bruggemann’s equation, that is 


DË = (e)X D; (4) 


where D; is the mass diffusivity of species i and x is the Brugge- 
mann coefficient of the porous layer. 


2.1.2. Electrochemical reaction (catalyst layers) 

In Eq. (3), the term S; stands for the source term of species i 
coming from the electrochemical reactions in the catalyst layers. 
Note that Hz and O2 are consumed by the reactions in the catalyst 
layers of anode and cathode, respectively, and H20 is the product 
of the reaction at the cathode. Therefore, the source terms for 
different species are 


Sm =- 5 

RA j- (5a) 
je 

So, = — 5b 

o = IF (5b) 

Sig a (Sc) 


where ja and je are the transfer current densities at the anode 
and the cathode, respectively, which are calculated based on the 
Butler—Volmer condition as 


. g Yp, aa Qaa F —Qa c F 
= x |ex ex 
Ja = Ja,ref Ya ref p RT n p RT n 


(6a) 


; F Yo, "e Qca F —Ac,cF 
= x |ex ex 
Je = Je,ref ocd p RT n P RT n 


(6b) 


where ja ref and je,ref are the reference transfer current densities, 
and ra and rc the concentration parameters, at the anode and the 
cathode, respectively. In addition, aq and a. are the anodic 
and the cathodic charge transfer coefficients for anode reaction, 
and æc,a, and Æc, the anodic and the cathodic charge transfer 
coefficients for cathode reaction. Yo, ref and Yu, ref are the ref- 
erence concentrations at which the exchange current densities 
were obtained. All the reference values of these coefficients are 


already given in Table 1. Note that the reference current den- 
sities, transfer coefficients, and reference concentrations have 
been set so as to fit the experimental data. 


2.1.3. Electronic conduction (carbon plates, GDLs, and 
catalyst layers) 

The carbon plate, the GDL, and the catalyst layers all serve 
as conductors for electric current. In the carbon plate and the 
GDL, there is no electrochemical reaction taking place. The 
electrochemical reactions only take place in the catalyst layers. 
Due to the electrochemical reactions activated in the catalyst 
layers, there exists a source term in the electronic conduction 
equation. Therefore, in these three components the electronic 
conduction equation may be expressed in general form as 


V - (ove) =¢ (7) 


where ¢ is the electronic potential and o°f denotes the effec- 
tive electronic conductivity. For the carbon plates, the effective 
electronic conductivity of is treated as a constant (op). For the 
porous gas diffusion and catalyst layers, the magnitude of the 
effective electronic conductivity is determined in terms of the 


porosity in accordance with the Bruggemann’s equation as 


ff 

ogpL = (1 — eapi)**- oGpL (8a) 
ff 

Oat = Eka OS, Cat (8a) 


It is noticed that the effective electronic conductivity of the cat- 
alyst layer (on) increases with the volume fraction of solid 
catalyst particles in catalyst layer (€s cat). In addition, the source 
term ¢ in Eq. (7) is zero for the carbon plates or the GDL with- 
out electrochemical reaction. When the catalyst layers is dealt 
with, ¢ is equal to the transfer current densities —j, and je for the 
anode and the cathode, respectively, which are calculated with 


Eqs. (6a) and (6b). 


2.1.4. Ionic conduction (catalyst layers and proton 
exchange membrane) 

The protons of hydrogen (H*) are generated in the anode 
catalyst layer and then migrate through the proton exchange 
membrane toward the cathode catalyst layer to react with oxygen 
and the returning electrons. Therefore, ionic conduction occurs 
in the catalyst layers and the proton exchange membrane. 

The ionic conduction equation for determination of the phase 
potential (g) in the catalyst layers is 


-V (EV Gcat) = 6 (9) 


where @car is the phase potential and I~ a is the effective ionic 
conductivity of the catalyst layers. Note that the source term ¢ is 
identical to that of Eq. (7). The effective ionic conductivity r 
is evaluated by using 


ff : 
Feat = EN Cat TN, Cat (10) 


It is important to note that since the content of Nafion in the 
catalyst layers help reduce the ionic resistance, the effective 
ionic conductivity of the catalyst layer cs) increases with the 
volume fraction of Nafion in catalyst layer (€N cat). 
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On the other hand, the ionic conduction equation for the 
membrane is expressed as 


V-UmV¢m) = 0 (1) 


where 9m is the phase potential and Im is the ionic conductivity 
of the membrane. The expression for ionic conductivity of the 
membrane has been proposed by Springer et al. [16] as: 


1 1 
Ty(T) = 1268 | —- — — 0.005239A — 0.00326 
m(T) ap | (= alk ) 
(12) 


where the water content of the membrane (A) is dependent on 
the dimensionless water activity in the membrane (a) accord- 
ing to experimental data. The activity is assumed uniform over 
the membrane, which is defined as the gas partial pressure of 
water at equilibrium with the membrane divided by the satura- 
tion pressure of water vapor at the operating temperature. That 
is 


Xy,0P 
a = —— 


(13) 
Psat 


where the saturation pressure of water vapor can be computed 
from 


logo Psat = —2.1794 + 0.02953(T — 273.15) — 9.1837 

x 107°(T — 273.15)? + 1.4454 

x 1077(T — 273.15) (14) 
with T in K. The correlation between À and a is derived as 


_ J 0.043 + 17.81a — 39.85a? + 36a? O<a<1 as 
~~) 144 14@= 5 l<a<3 


Based on the solutions for the phase potential (g), the distribu- 
tion of local ionic current density (i;) may be further determined 
with i = —I'V@g. 


2.1.5. Nafion loading (catalyst) 

As mentioned earlier, Nafion is an electronic insulator. Hence, 
when it is introduced into the catalyst layer, it decreases the elec- 
tronic conductivity of the catalyst layer. However, on the other 
hand, since it is a good conductor for the hydrogen protons, 
Nafion loading can elevate the ionic conductivity of the catalyst 
layers. Besides, the introduction of Nafion into the catalyst lay- 
ers changes simultaneously the porosity and the gas permeability 
of the layers so that it affects the transport of the gas species. 
As Nafion is loaded, total volume of the catalyst layer (v;) may 
be divided into three major portions, namely, the volume occu- 
pied by the solid catalyst particles (vs), the volume occupied by 
Nafion (vy), and the volume of the void space (vy). The volume 
fractions of the three portions are denoted by &s cat, €N,Cat, and 
EV,Cat, respectively, and are defined as 


US 


ES,Cat = — (16a) 
Ut 
v 

EN,Ca = ~ (16b) 
Ut 


and 


Ev,Cat = l — EC,Cat — EN,Cat (16c) 


It is noted that the present model does not include the effects 
of heat and liquid water accumulation. A more powerful model 
will absolutely be desired when necessary. 


2.1.6. Solution domain and boundary conditions 

For the pattern of straight parallel channels shown in Fig. 1, as 
the number of gas channels is large enough, a periodic transport 
phenomenon is expected to develop from channel to channel. 
Thus, one only needs to deal with the solution domain of a 
single module indicated by the dashed lines. In this study, the 
width of a module (/p) is set to be 1 mm. Therefore, for the faces 
of solution domain at x=0 and 1 mm, a symmetric boundary 
condition for U , Y, ġ, and ọ is prescribed. 

Solid walls of the gas channels or the carbon plates are imper- 
meable, and hence, on the solid walls, 9Y;/ð3n=0 and Ù =0. 
In addition, since both the gas channels and the membrane 
are insulators against electrons; at the boundaries of these two 
components, the normal gradient of the electric potential (@) is 
vanished. Furthermore, the ionic conduction takes place only 
inside the catalyst layer and proton exchange membrane, and 
the protons are not allowed to transfer into other components; 
therefore, on the faces of the catalyst layer and the membrane 
the normal gradient of phase potential is assigned to be zero. 

As to the inlet conditions, the inlet velocities for the anode gas 
(hydrogen) and the cathode gas (air) in the gas channel, V and 
Vo, are fixed at 0.3 and 0.5 m s7}, respectively. The inlet mass 
fractions for all the gas species are specified, and their values 
are provided in Table 1. Note that both the fuel and the oxidant 
gases enter the fuel cell at 100% relative humidity. On the other 
hand, for the boundary conditions at the exit, it is assumed that 
the flow becomes fully developed downstream. Thus, the exit 
boundary conditions may be described as 

dv OY; ə dy 


“=0= ay ay By ay 


=0 (17) 


2.2. Numerical methods 


These above set of conservation equations along with the 
property equations are solved by adopting a finite-volume 
scheme on structured grids within the framework of the com- 
mercial computational fluid dynamics code (CFD-ACE+). A 
self-developed subroutine has been built for calculating the 
values of the properties of the catalyst layers, including the effec- 
tive electronic conductivity, the effective ionic conductivity, and 
effective mass diffusivity of gas species (f, r, and D) for 
various Nafion loading conditions. The Nafion loading condition 
is specified by adjusting relatively the values of volume fractions 
of Nafion loaded, solid catalyst particles, and void space inside 
the catalyst layers (E€N,Cat, ES,Cat, and £vy,Cat). The self-developed 
subroutine is incorporated with the computational fluid dynam- 
ics code through a Python-interface connection. The message 
of any change in the Nafion content is firstly used by the self- 
developed subroutine to update the values of the properties of 
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Table 3 
Specification of single cell module and test conditions in experiments 


22cm x 22cm 

14cm x 14cm 
Five-layer (Dupont) 
Composite graphite plate 
(3-mm thickness) 
T=80°C and RH= 100% 
T=80°C and RH= 100% 
1500 sccm 

2500 sccm 


Cross-sectional area 
Active area of membrane 
Membrane 

Carbon plates 


Inlet conditions of hydrogen 
Inlet conditions of air 
Volumetric flow rate of hydrogen 
Volumetric flow rate of air 


the catalyst layer, and then the updated values of the catalyst 
layer properties can be transmitted automatically to the compu- 
tational fluid dynamics code through the interface in order for 
solution of the three-dimensional transport phenomena taking 
place in the PEMFC. When necessary, the CFD code solutions 
may also be sent back to the self-developed subroutine through 
the interface for further computation. 

In an average case, the structured grids containing approx- 
imately 100,000 to 200,000 control volumes are used. A 
grid-independence check has been performed to ensure the 
accuracy of the numerical solutions. At each iteration step, the 
iteration for the solutions for U , Y, ġ, and ¢ is progressed until 
the following convergence criterion is satisfied: 

CAA mie <10 (18) 
yo 
i,j,k 
where W represents each of the variables U , Y, @, and g; l the 
iteration number; and i, j, and k represent the spatial indices for 
the x-, y-, and z-directions, respectively. 

The average current density (I) of fuel cells is calculated by 
integration of local electric current density distribution on the 
carbon plate surface at z=0 divided by the carbon plate surface 
area as 


1 A 
=3 f is- dA (19) 


where i, is the local electric current density on the carbon plate 
surface and A is the carbon plate surface area. In the mean time, 
the gas channel width fraction (A) is defined by the ratio of the 
gas channel width (/c) to the width of a module (/p), that is 


=s (20) 


3. Experiments with single-cell PEMFC 


A single-cell fuel cell module with cross-sectional area 
of 22cm x 22cm and active area of 14cm x 14cm has been 
installed and tested to partly verify the numerical predictions 
yielded by the commercial computational fluid dynamics code. 
The specification of the single cell module and the test conditions 
are given in Table 3. The purpose of the experiments is simply 
to test the validity of the CFD code; therefore, Nafion is not yet 
introduced into the catalyst layers of the MEA. The operating 
and geometrical conditions in simulation are consistent with the 


Fig. 2. Experimental apparatus (left) and single-cell fuel cell module (right). 


experimental module. Two composite graphite plates of 3-mm 
thickness are used as the flow fields and current conductors. A 
five-layer membrane exchange assembly (MEA) provided by 
Dupont Inc. is placed at the center between the two composite 
graphite plates. The membrane used in this MEA is Nafion® 
112 and the gas diffusion layers are made of carbon fiber paper. 
Before entering the channels, the fuel and oxidant gases are 
humidified to maintain at 100% relative humidity and 80°C 
temperature. 

In addition, to obtain the polarization curves with the fuel cell, 
a performance test system has also been setup. In the exper- 
iments, the volumetric flow rates of the hydrogen and the air 
are typically maintained at 1500 and 2500 sccm, corresponding 
to V,=0.3ms~! and V; = 0.5 m s7!, respectively. The fuel cell 
temperature is maintained at 343 K. Photographs of the experi- 
mental apparatus and the single-cell fuel cell module are shown 
in Fig. 2. 
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Fig. 3. Polarization curves: comparison between predictions and experiments. 
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Fig. 3 shows a comparison between the numerical predic- 
tions and the experiment data for the polarization curves of the 
single-cell fuel cell module without Nafion loading. The oper- 
ating conditions have been described in the preceding section. 
It is observed that the numerical predictions closely agree with 
the experiments, especially as 1<3800 A m~?. The difference 
between the two sets of data becomes appreciable only in the 
high-current-density regime. In this regime, the experimental 
data appear to experience a greater drop in the operating voltage 
than the predicted. The greater voltage drop is attributed to the 
concentration loss due to insufficient supply of O2 and H? in 
experiments at the cathode and the anode, respectively which 
is remarkable in the high current density situations. In Fig. 3, 
it is found that at V=0.7 V the numerical and the experimental 
data for the current density are in close agreement (at approx- 
imately 2500 Am“). Thus, the typical value of the operating 


x=Imm 
(a) (b) 


x=0 mm 


Cathode 


voltage of the fuel cell is fixed at 0.7 V in the present study to 
evaluate the Nafion loading effects. Note that in the following 
parametric study, the operating temperature is assigned to be 
323 K. However, the experiments have been performed at dif- 
ferent temperatures. The purpose of the experiments is simply 
to test the validity of the CFD code; therefore, only the data at 
343 K are shown in Fig. 3. 


4. Results and discussion 


Parametric study of the effects of the parameters listed in 
Tables | and 2 is now carried out to evaluate the performance of 
PEMFC under various combinations of these parameters. 

Nafion loading effects on the distributions of mass frac- 
tions of hydrogen and oxygen at the anode and the cathode, 
respectively, are shown in Fig. 4. Numerical solutions on 


(c) 


Fig. 4. Mass fraction distributions on cross-section at y= 25 mm with different volume fractions of Nafion in catalyst layer, for V=0.7 V, egpL =0.5, €s,cat = 0.885, 
lc =5 x 1074 m, tapL =3 x 1074 m, tm = 1.76 x 10-4 m, h=1 x 1073 m, and A=0.5. (a) en,cat = 0.1; (b) EN.cat = 0.2; (C) EN,Cat = 0.3. 
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Fig. 5. Planar distribution of ionic current density at the center plane of the cathode catalyst layer with different volume fractions of Nafion in catalyst layer, 
for V=0.7 V, eapL =0.5, £s Cat =0.885, Ic =5 x 1074m, tapL =3 x 1074m, fm =1.76 x 1074m, h=1 x 107°? m, and A=0.5. (a) EN,cat =0.1; (b) £N,cat =0.2; (c) 


EN,Cat = 0.3. 


cross section at y=25mm with different volume fractions of 
Nafion in catalyst layer (€N,Ca=0.1, 0.2, and 0.3) are plot- 
ted, for V =0.7 V, egpt = 0.5, £s,Cat = 0.885, tGDL =5 x 1074 m, 
fm = 1.76 x 107-4 m, h=1 x 107? m, and A =0.5. It is expected 
that the introduction of Nafion reduces the porosity of the cat- 
alyst layers; therefore, the mass fractions of the reactant gases 
in the catalyst layers are obviously decreased by elevating the 
volume fraction of Nafion. However, since Nafion loading can 
elevate the ionic conductivity of the catalyst layers, the current 
density may be increased by the introduction of Nafion. For 
example, Fig. 5 displays the planar distribution of the ionic cur- 
rent density at the center plane of the cathode catalyst layer in the 
same conditions of Fig. 4. It is clearly seen that the uniformity 
of the ionic current density distribution on the plane is signifi- 
cantly improved by increasing the volume fraction of Nafion in 
catalyst layer. In accordance with the discussion concerning the 


gases transport and the ionic conduction discussed herein, an 
optimal Nafion content in the catalyst layer is indeed expected 
and desired. 

Given the operating conditions of V=0.7 V, egpL =0.5, 
Es,Cat =0.588, tcpL=5 x 1074m, fm = 1.76 x 1074m, h=1 x 
107° m, and A =0.5, one may examine the dependence of power 
density of the PEMFC (P) on the volume fraction of Nafion in 
the catalyst layer (€n cat). The predicted results are plotted in 
Fig. 6. It is clearly seen that the magnitude of the power density 
first increases with £N,Cat to reach its optimal value of approx- 
imately 1812 W m~” at EN,Cat = 0.3, and then it descends when 
the value of én cat over the optimal one. 

The data shown in Fig. 6 are collected as the volume fraction 
of solid catalyst particles in the catalyst layer (€s cat) is fixed 
at 0.588. On the other hand, what plotted in Fig. 7 is the infor- 
mation concerning the dependence of power density on én cat 
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Fig. 6. Power density as a function of volume fraction of Nafion in catalyst layer 
at V=0.7V, ecpL =0.5, £s,ca =0.588, Ic=5 x 10-4 m, top. =3 x 1074m, 
tm= 1.76 x 1074m, h=1 x 107° m, and A=0.5. 


at a constant volume fraction of void space in the catalyst layer 
(€v,cat). Note that the solid volume fraction is reduced as more 
Nafion is added while the void space stays constant. Again, the 
curve exhibits a peak value as £N,Cat is varied from 0.1 to 0.8. 
Results proved in Fig. 7 illustrates that in this case the maxi- 
mum power density can be yielded when the volume fraction of 
Nafion in the catalyst layer is maintained at 0.45. 

Fig. 8 shows the variation of power density with the vol- 
ume fraction of void space in the catalyst layer for V =0.7 V, 
€GDL=9.5, EN,Cat =0.3, Ic =5 x 10-4 m, GDL =3 Xx 10-4 m, 
tm = 1.76 x 10-4 m, h=1 x 10-7 m, and A=0.5. The volume 
fraction of void space in the catalyst layer (€y,cat) is also referred 
to as the porosity of the catalyst layer. For this case, the ionomer 
content stays constant and in order to increase porosity only the 
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Fig. 7. Power density as a function of volume fraction of Nafion in catalyst layer 
at V=0.7V, egpL =0.5, evicar=0.112, Ic=5 x 10-4 m, top. =3 x 1074m, 
tm= 1.76 x 10-4m, h=1 x 107° m, and A=0.5. 
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Fig. 8. Power density as a function of volume fraction of void space 
in catalyst layer at V=0.7V, egpL =0.5, £s Cat =0.588, Ic=5 x 1074 m, 
tapL =3 x 1074 m, fm = 1.76 x 10-4 m, h=1 x 107° m, and A=0.5. 


solid phase is decreased. An increase in £y,Cat is favorable for 
sufficient gas supply; however, it reduces the electronic conduc- 
tion. For the case considered in this figure, it is observed that the 
optimal power density is reached at ey,cat =0.125. 

Fig. 9 further conveys the effects of gas channel width 
fraction (A) on the performance of the fuel cell at V=0.7 V, 
eGpL=0.5, €v.cat=0.112, en.car=0.3, top. =3 x 1074m, 
tm=1.76 x 10-4m, and h=1 x 107°m. It is observed that 
indeed the channel width fraction has subtle influence on the 
cell performance. Larger channel width fraction means a larger 
area occupied by the gas channels and more reactant gases are 
supplied; however, it also means thinner ribs and hence smaller 
area for electron conduction. Therefore, the ohmic impedance 
increases with the channel width fraction. With these two oppos- 
ing effects, there must exists an optimal value of the gas channel 
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Fig. 9. Power density as a function of channel width fraction at V=0.7 V, 
EGDL = 0.5, €v.cat = 0.112, EN,Cat = 0.3, tap =3 x 1074 m, tm = 1.76 x 1074 m, 
andh=1 x 1073 m. 
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Fig. 10. Power density as a function of thickness of gas diffusion 
layer at V=0.7V, egpL=0.5, vat =0.112, encat=0.3, Ic =5 x 1074 m, 
tm = 1.76 x 1074 m, h=1 x 107° m, and A=0.5. 


width fraction. In Fig. 9, the results show that the largest 
power of 1816 W m7? is obtained at the optimal gas channel 
width fraction of A=0.4 under the conditions considered in 
this case. 

Plotted in Fig. 10 is the variation in the power density with the 
thickness of gas diffusion layer (tgp) for V =0.7 V, egpL = 0.5, 
€v.Cat = 0.112, en,cat = 0.3, Ic =5 x 1074 m, tm = 1.76 x 1074 m, 
h=1x10-*m, and A=0.5. A thicker gas diffusion layer 
may allow the gas to disperse uniformly in the entire gas 
diffusion layer; however, it also tends to increase the elec- 
tric resistance against the conduction of electrons. It is seen 
that the optimal performance of the fuel cell is 1828 W m~? 
at tcpL =200 um. In addition, the porosity of the gas diffu- 
sion layer (gp) is also one of the major factors affecting 
the fuel cell performance. Similarly, an increase in €GpL is 
favorable for sufficient gas supply; however, it also tends 
to increase the electric resistance. For the case at V=0.7 V, 
Ev,Cat = 0.112, en,cat = 0.3, lc =5 x 1074 m, top, =3 x 1074 m, 
tm = 1.76 x 10-4m, h=1x 10-3m, and A=0.5, the results 
shown in Fig. 11 illustrate that the peak value of the power 
density can be achieved at €gpL = 0.3. 

Attention is now drawn to the effects of the thickness of 
the proton exchange membrane on the performance of fuel 
cell. Fig. 12 displays the dependence of power density on 
the thickness of the membrane (tm) at V=0.7 V, egpr = 0.5, 
€v.Cat = 0.112, en,cat = 0.3, lc =5 x 1074 m, top, =3 x 1074 m, 
h=1 x 107°? m, and A=0.5. It is noticed that within the range 
between tm=50 and 500 um, the power density of the fuel 
cell is monotonically decreased by increasing the thickness of 
the membrane, and hence, no optimal membrane thickness is 
observed. This reflects the fact that a thicker membrane possess- 
ing higher mechanical strength is also a membrane possessing 
higher ionic resistance. However, the other important aspect of 
the membrane thickness effects is that the cross-over current, 
which reduces the performance of a fuel cell, increases with 
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Fig. 11. Power density as a function of porosity of gas diffusion layer 
at V=0.7V, evcat=0.112, encar=0.3, lc=5 x 1074m, top. =3 x 1074m, 
tm = 1.76 x 10-4m, h=1 x 107° m, and A=0.5. 


decreasing membrane thickness. Therefore, it is not true that 
thinner membranes are always better. 

It is important to mention that in order to be able to report 
optimal conditions for more than one design variable, a mul- 
tivariable numerical optimization problem (e.g. Lin et al. [10], 
Song et al. [13,14], and Secanell et al. [15]) is necessary. A good 
example of the need for such an approach is the difference in 
optimum Nafion content between Figs. 6 and 7 for which a sin- 
gle parameter was changed. Further, the optimum channel width 
fraction is 0.4 only for a GDL with 0.5 porosity and a catalyst 
layer with 0.112 porosity and 0.3 volume fraction of ionomer. 
This value might be completely different for a catalyst layer with 
0.8 porosity. 
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Fig. 12. Power density as a function of thickness of the membrane at V =0.7 V, 
ecpL=0.5, €vcar=0.112, en,ca=0.3, Ic=5x 1074m, top, =3 x 1074m, 
h=1 x 1073 m, and A=0.5. 
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5. Conclusions 


This study is concerned with the Nafion loading effects 
by numerical simulation using a full-cell three-dimensional 
numerical model. In addition, change in the performance of 
the PEMFC accompanying various combinations of geometric 
parameters is also investigated. A 22 cm x 22 cm single-cell fuel 
cell module with 14cm x 14cm active area has been installed 
and tested. Part of the numerical results is compared with exper- 
iments for validation. It has been observed that the numerical 
predictions closely agree with the experiments, especially as 
1<3800A m7. 

A self-developed subroutine has been built for calculating 
the values of the properties of the catalyst layers, including the 
effective electronic conductivity, the effective ionic conductiv- 
ity, and effective mass diffusivity of gas species for various 
Nafion loading conditions. The self-developed subroutine is 
incorporated with a computational fluid dynamics code through 
a Python-interface connection to deal with the three-dimensional 
mass, momentum, and species transport phenomena as well as 
the electron- and proton-transfer process taking place in the 
PEMFC. 

For the particular cases considered in this study, the optimal 
conditions of the volume fraction of Nafion in the catalyst layer 
(EN,Cat), the volume fraction of void space in the catalyst layer 
(Ev,Cat), the gas channel width (/c), and the thickness and poros- 
ity of the gas diffusion layer (tGpL and eGpL), which lead to a 
higher fuel cell power density, are investigated. Power density of 
the fuel cell is monotonically decreased by increasing the thick- 
ness of the membrane (tm), and hence, no optimal membrane 
thickness is observed. Note that the suggested optimal param- 
eters are obtained under the base case conditions considered in 
the present study. When the conditions are changed, the opti- 
mal values of these parameters may be altered. Therefore, in 
order to be able to report optimal conditions for more than one 


design variable, a multivariable numerical optimization problem 
is necessary. 
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